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Abstract. Two questions that naturally arise in A-body 
simulations of stellar systems are: 

1. How can we compare experiments that employ differ- 
ent types of softened gravity? 

2. Given a particular type of softened gravity, which 
choices of the softening length optimize the faithful- 
ness of the experiments to the Newtonian dynamics? 

We devise a method for exploring the dynamical effects of 
softening, which provides detailed answers in the case of 
2-D simulations of disc galaxies and also solves important 
aspects of the 3-D problem. In the present paper we focus 
on two applications that reveal the dynamical differences 
between the most representative types of softened gravity, 
including certain anisotropic alternatives. Our method is 
potentially important not only for testing but also for de- 
veloping new ideas about softening. Indeed, it opens a 
direct route to the discovery of optimal types of softened 
gravity for given dynamical requirements, and thus to the 
accomplishment of a physically consistent modelling. 

Key words: gravitation - methods: analytical - methods: 
numerical - galaxies: evolution - galaxies: kinematics and 
dynamics - galaxies: spiral 



1. Introduction 

In ./V-body simulations of stellar systems the gravitational 
interaction is modified for curing the Newtonian diver- 
gence at short distances. Basically, such modifications are 
introduced through a soft cut-off: the softening length s. 
However, the precise form in which they are implemented 
can vary. Since gravity plays a fundamental role in these 
systems and the gravitational interaction is modified pre- 
cisely where it becomes singular, the dynamical effects 



of softening should be well understood when designing 
experiments and interpreting their results. This dynam- 
ical problem has recently stimulated considerable inter- 
est (e.g., Hernquist & Barnes 1990; Hernquist & Ostriker 
1992; Kandrup et al. 1992; Pfenniger 1993; Pfenniger & 
Friedli 1993; Gurzadyan & Pfenniger 1994; Romeo 1994, 
hereafter Paper I; Byrd 1995; Gerber 1996; Merritt 1996; 
Weinberg 1996; Sommer-Larsen et al. 1997; Theis 1997; 
see also Goodman et al. 1993; Farouki & Salpeter 1994). 
For extensive overviews see the above-mentioned Pfen- 
niger & Friedli (1993), Gurzadyan & Pfenniger (1994) and 
Paper I. 

In Paper I we have investigated the stability problem 
in the case of 2-D models with Plummer softening, which 
are commonly employed in simulations of disc galaxies. 
The basic message is that the effect of softening becomes 
strongly artificial for s ^ X/2ir, A being the typical ra- 
dial wavelength, which means half an order of magnitude 
below the expected value. The major results are summa- 
rized in the form of a criterion of approximate physical 
consistency for s and a stability criterion for the Toomre 
parameter. (Other important aspects of the stability prob- 
lem have been considered by Byrd 1995.) 

In the present paper we carry out five extensions, as is 
discussed below. 

1. We generalize the stability analysis of Paper I to an 
arbitrary isotropic form of softening. This is a natural 
extension since types of softened gravity different from 
the standard Plummer softening are becoming more 
and more commonly employed (e.g., Combes et al. 
1990; Parous et al. 1993; Shlosman & Noguchi 1993). 
In particular, the alternatives proposed by Hernquist 
& Katz (1989) and Pfenniger & Friedli (1993) reflect 
an interesting idea, viz. that softening should be as lo- 
calized as possible since there is no clear reason for 
modifying the gravitational interaction at long dis- 
tances, and it is tempting to explore its dynamical con- 
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sequences. From a more general point of view, this and 
the following extensions provide the tools for compar- 
ing experiments that employ different types of softened 
gravity. 

2. We investigate the implications of our stability analy- 
sis for the classical relaxation problem (Rybicki 1972; 
White 1988). Relaxation and stability 1 are intimately 
related in self-gravitating systems, and even simple 
treatments reveal their strong coupling through ran- 
dom motions. On the other hand, the contribution of 
velocity dispersion to the relaxation time has been un- 
derstood only in part and, because of that, the classi- 
cal argument favouring the choice of large values of s 
is wrong. We revise this argument and conclude that 
neither small nor large values of s are convenient. Sur- 
prisingly, there exists an intermediate choice of s that 
optimizes the 'dynamical resolution' of the model, i.e. 
its faithfulness in simulating the dynamics of 3-D discs 
with Newtonian gravity, especially in situations near to 
the stability threshold. We identify the optimal char- 
acteristics, and show how to evaluate them for a given 
type of softened gravity. In addition to investigating 
this aspect of the relaxation problem, we explain how 
effectively softening reduces noise on various scales. 

3. We complete the examination of 2-D models with 
isotropic softening by investigating the equilibrium 
problem for an axisymmctric state with epicyclic mo- 
tions. In particular, we explain how significantly the 
circular speed and related quantities deviate from their 
Newtonian behaviours at various distances from the 
centre. 

4. We consider 3-D models with isotropic softening and 
examine two limiting cases: discs and the simple, yet 
instructive, Jeans problem. An extension to 3D has 
been encouraged by Friedli (1994) and Junqueira & 
Combes (1996, see the interesting remarks in Sect. 
2.2). Real stellar systems have several gravitationally 
interacting components. Both iV-body simulations and 
theoretical works are forced to use simplified models, 
which do not necessarily provide faithful representa- 
tions of the complexity of such systems. Our motiva- 
tion is to understand the basic differences between the 
dynamical effects of softening in 3D and 2D, in the 
presence of a single stellar component (3-D vs. 2-D 
discs and Jeans problem vs. discs). In particular, we 
point up the strong modifications introduced by a ho- 
mogeneous geometry and the absence of rotation. 

5. We complete the examination of 3-D models by 
discussing the basic dynamical effects of softening 
anisotropy. The idea underlying this extension is that 
softening should be anisotropic in simulations of stellar 



1 In this and similar contexts, 'stability' should be under- 
stood in the general sense of 'stability properties'; it is not 
implied that the system is stable. The same applies to the use 
of 'relaxation' and 'equilibrium'. 



systems where significantly higher spatial resolution 
is required along a certain direction, such as in disc 
galaxies. The alternative family of softening recently 
proposed by Pfenniger & Friedli (1993) reflects such an 
important idea, and it is tempting to explore the dy- 
namical relations between its members. (An analogous 
idea has been discussed in the context of smoothed 
particle hydrodynamics by Shapiro et al. 1994, 1996 
and Fulbright et al. 1995; Hcrnquist, private commu- 
nication, has remarked that in models with anisotropic 
smoothing there may be a significant tendency for an- 
gular momentum not to be conserved.) 

These extensions all together form a method for exploring 
the dynamical effects of softening in A-body simulations 
of stellar systems. Our method is described in Sect. 2, and 
is structured as in the previous discussion. The two appli- 
cations mentioned in the same context are shown in Sect. 
3 (see also Appendix A). The conclusions of this paper 
are drawn in Sect. 4, where we present our contribution 
in a more general perspective and motivate future appli- 
cations. 

2. Method 

2.1. 2-D models with isotropic softening 

2.1.1. Stability 

In 2-D discs with isotropic softened gravity, a given 
surface-density perturbation £i(i?, t) induces a potential 
perturbation 



$i(R,t) = -G fip s (\R - R'\) S 1 (i?',t)d 2 J R' , 



(1) 



where —Gni(p s (R) is the point-mass potential, which de- 
pends on the softening length s [ip s (R) is also known as 
the softening kernel]. In analysing the Poisson equation 
(1), we adopt the lowest-order WKBJ approach, as in Pa- 
per I, but with a different albeit asymptotically equivalent 
spectral representation. In simple terms, the new feature 
consists in considering perturbations with a radial depen- 
dence gi(R) = giJo(kR) rather than gi(R) = g\e lkR , 
where k is the radial wavenumbcr and J„ denotes the 
Bessel function of the first kind and order v (see, e.g., 
Abramowitz & Stcgun 1972). The Bessel-Hankel represen- 
tation is more convenient than the Fourier representation 
because it allows factorizing the convolution in the Pois- 
son equation directly, without any assumption on the form 
of <p a (R): 



*i 



-2nG H [cp s (R)}(k) Ex, 



(2) 



where 7i v denotes the Hankel transform of order v (see, 
e.g., Sneddon 1972; Bracewell 1986): 



/■OO 

H v [g{R)]{k)= I g(R)J v (kR)RdR. 
Jo 



(3) 
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The reduction factor. Equation (2) admits of a simple 
interpretation: the potential perturbation induced by a 
given surface-density perturbation is weakened by a fac- 
tor 



S 



-Ho[MR)Kk) 
H [MR)](k) 



(4) 



—Grrup-^(R) being the Newtonian point-mass potential. 
Correspondingly the contribution of self-gravity to the 
dispersion relation is weakened through a reduction of the 
active unperturbed surface density by the same factor. 
Thus, S provides complete information about the effect 
of isotropic softening on the stability of 2-D discs, and 
its evaluation is the starting-point of our method. There 
are two complementary routes to S. One is to calculate it 
analytically from: 



S(\k\s) = |*| H [vs{R)](k) 



(5) 



(use, e.g., the comprehensive table of integrals by Grad- 
shteyn & Ryzhik 1994). This equation follows directly 
fromEq. (4) , being ip N (R) = l/i?and H [l/R]{k) = l/\k\. 
The other is to compute it numerically from: 



S(\k\s) = l-H 1 



ff-f'W 



(1*1)- 



(6) 



—Gm 2 f s (R) being the point-mass force (use, e.g., the ex- 
tensive and well-documented NAG library; for algorithms 
of fast Hankel transform see, e.g., Gueron 1994; van Veld- 
huizen ct al. 1994). This equation is more convenient than 
the original Eq. (5) because the slow decay of the oscilla- 
tory integrand as R — ► oo is speeded up. Such an improve- 
ment is obtained with simple tricks, viz. integrating by 
parts and singling out the Newtonian behaviour of f s (R) 
at long distances. 

The effective scale height. Once <S(|fc|s) has been evalu- 
ated, the second step is to analyse its behaviour at small 
\k\s since stability is basically a large-scale property. In or- 
der to understand the general features of this behaviour, 
we can equivalently start from Eq. (5) or (6) and use 
techniques of asymptotic expansion of integrals (see, e.g., 
Bender & Orszag 1978). In Eq. (5), we should first sin- 
gle out the Newtonian 1/R -dependence of (fi s (R) at long 
distances, and then expand Jo(fci?). In Eq. (6), we can di- 
rectly expand Ji(|fc|_R). The result is that S ~ 1 — \k\h, 
where in each case 



[l-Rip s (R)]dR 



[l-R 2 f s (R)]dR. (7) 



This quantity is positive in types of softened gravity of 
practical interest, and has an important dynamical mean- 
ing. In fact, a comparison with the reduction factor of 
3-D discs with Newtonian gravity (Shu 1968; Vandervoort 
1970; Romeo 1992) shows that softening mimics thickness 



on large scales and h has the effect of a scale height, as 
far as density waves are concerned (for bending waves see 
Masset & Tagger 1996). As an alternative to Eq. (7), h 
can be evaluated from the conversion factor 



- = -S'(\k\s = 0). 

s 



(8) 



Characteristics. The third and last step is to extract 
detailed information concerning the stability properties. 
This part of the method has been described in Paper I 
and can be generalized without special difficulties. So in 
the following discussion we briefly introduce the basic con- 
cepts and point out the major results of the stability anal- 
ysis. We adopt the basic fluid description with the follow- 
ing scaling and parametrization: 



A = - — , where At 
Ax 



2tt _ 4tt 2 G£ 

fcx K 2 



Q = (Safronov-Toomre parameter) . 

rj = sfcx (softening parameter) , 

£ = /ifcx (effective thickness parameter) , 



(9) 

(10) 

(11) 
(12) 



where A is the radial wavelength of the perturbation, Ax 
is the Toomre wavelength, £ is the unperturbed surface 
density, k is the epicyclic frequency and c is the planar 
sound speed. The stability properties are described by the 
marginal stability curve, i.e. the dispersion relation for 
marginally stable perturbations viewed as a curve in the 
(A, Q 2 ) plane for a given r\: 



Q 2 + 4A [A - S(r]/X)] =0. 



(13) 



In particular, the stability level is measured by the effec- 
tive parameter 



Qe 



Q 



(14) 



where the threshold Q{rf) corresponds to the square root 
of the global maximum of the marginal curve, and the typ- 
ical radial wavelength corresponds to the location A max (77) 
of this maximum. The values of 77 and the related quanti- 
ties that characterize the stability properties are presented 
below. 

• The conversion factor C,/ri — —S'(r)/\ = 0), now ex- 
pressed in dimensionless form, has already been dis- 
cussed (cf. the effective scale height). 

• The safety threshold for approximate physical consis- 
tency is ?7 sa f = i (C/r/y 1 . For 77 <, r? sa f e , i.e. C ^ 5, 
softening mimics the effect of thickness. In particular, 
Q»l-2( and A max ~ \. For 77 ^ ?7 sa f , softening 
causes artificial stabilization and a moderate degree of 
'blueshift'. 



Alessandro B. Romeo: Dynamical effects of softening in iV-body simulations of disc galaxies 



• The critical value ?7 cr it is such that Q(?7crit) = 0. For 
r] > ?7 cr jt , the fundamental meaning of velocity disper- 
sion becomes obscure because the stability level is no 
longer actively controlled by Q c g. This fact has seri- 
ous dynamical implications: it artificially precludes the 
possibility of simulating regimes of normal spiral struc- 
ture, which require fine-tuned choices of the stability 
level. 

• The critical radial wavelength is A cr ;t = A max (??crit)- Its 
deviation from A max (0) = i indicates the sensitivity of 

Amax(?7) for Tj < Tfcrif 

2.1.2. Relaxation vs. stability 

Let us investigate the implications of our stability analysis 
for the classical relaxation problem (Rybicki 1972; White 
1988; see also Hockney & Eastwood 1988). We first sketch 
the basic ideas. The Rybicki- White relaxation time £rw is 
proportional to the softening length s, the cube of velocity 
dispersion a and the number of computer particles N. In 
turn, a is proportional to the Safronov-Toomre parameter 
Q. So it seems that large values of s are convenient, since 
they result in a long t-gw for given Q and N. The idea 
underlying this argument is that a given Q corresponds 
to a constant stability level as s varies. But it is not so. 
In fact, the stability threshold Q depends on s and the 
level is measured by the effective parameter Q e g = Q/Q. 
Thus large values of s are not convenient at all, because 
they correspond to a low Q and result in a short £rw fo r 
given Q ff and N. Furthermore, if %w is short for both 
small and large values of s, there must exist an interme- 
diate choice of s that maximizes few ■ Does it satisfy the 
criterion of approximate physical consistency? The answer 
is: Yes, it does. The identification of the optimal charac- 
teristics of 2-D discs with isotropic softened gravity is the 
focal point of our method. 

The corrected Rybicki- White relaxation time. Originally, 
few was derived by adopting a simple two-body treat- 
ment, and assuming that f s (R) is of the type: f s (R) = 
for R < s, f s (R) — 1/R 2 for R > s. A generalization of 
few to an arbitrary isotropic f s {R) is easy to derive and 
useful for comparing the effects of different types of soft- 
ened gravity. The resulting t RW is equal to the original 
one multiplied by a correction factor C, which is one of 
the relaxation characteristics. Referring to the introduc- 
tory discussion, we express few in a form that splits the 
various contributions: 



where r(r]) measures the 'relaxation level' as 77 varies, M^ 
is the disc mass and b is the impact parameter. In addi- 
tion, the original unspecified a is interpreted as the radial 
velocity dispersion c (different specifications would only 
modify the proportionality factor). 



Fig. 1. The relaxation level t(ji) — nQ (n) of 2-D discs with 
Plummer softened gravity, where n is the softening parameter 
and Q{n) is the stability threshold. The contribution of Q to 
the ^-dependence of r has major dynamical implications, as is 
pointed up in Sect. 2.1.2 



Optimal characteristics. We now identify the values of rj 
and the related quantities that optimize the relaxation 
and stability properties. Figure 1 shows the behaviour of 
t(ji) in the case of Plummer softening, but the following 
discussion is general. 

• The optimal relaxation level r op and the optimal choice 
rj p are such that r op = t(t7 op ) = max{r(?7)}. The sim- 
ple analytical approximation r) op w g (C/v)^ 1 shows 
that 7/ p < 77 safc . 

• In 3-D discs with Newtonian gravity, the temperature 
anisotropy corresponding to the stability threshold for 
the optimal value Cop ~ g would be <5 op ~ 5, which 
means a realistic vertical-to-radial velocity dispersion 
ratio c z /c ps 0.6 [cf. Paper I, Eq. (16) and Fig. 3]. 

• The values r^pi and rj op 2 such that, e.g., r(?7opi) = 



l"(??op2) = 

of »?, ?7o P i 



r op specify a range of convenient choices 
V SS ? 7op2) and complete the information 



few = C ■ r{ri) Q 3 cS N 



2tt 3 G 2 £ 3 
n 5 M d 



(15) 



concerning the optimal characteristics. The suggested 
definition is natural and especially meaningful because 
it turns out that r) op 2 ~ ??safe- 



Cs 



[bfs(b)} 2 db 



r(v)=vQ 3 (v)- 



(16) 



What about the reduction factor? A delicate aspect of the 
relaxation problem that has not been considered in the 
previous discussion concerns the effects of collective inter- 
(17) actions between particles and self-consistent fluctuations 
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on the dynamical evolution of the system (e.g., Romeo 
1990 and references therein; Weinberg 1993; Zhang 1996). 
A thorough treatment of collective effects would demand 
titanic efforts even in simpler models (cf. Weinberg 1993). 
On the other hand, useful information is already contained 
in <S(|fc|s) since the diffusion properties are determined 
by the dispersion relation. In particular, the behaviour of 
<S(|fc|s) at large \k\s (S <C 1) shows how effectively soften- 
ing suppresses small-scale fluctuations, which represent an 
important source of noise in 2-D models (cf. White 1988; 
Schrocder & Comins 1989). 

2.1.3. Equilibrium 

The reduction factor. The equilibrium problem for an ax- 
isymmetric state with epicyclic motions can be solved by 
using the technique of Hankcl transforms, which has al- 
ready been introduced in the stability analysis. The basic 
idea is to Ho-transform the Poisson equation twice: once 
for factorizing the convolution of <p s (R) and X(-R), and the 
other time for recovering $(-/?) [or, in the inverse problem, 
S(i?)] from its transform. This is a natural generalization 
of the approach adopted by Toomre (1963) in the case with 
Newtonian gravity, where the Poisson equation can be ex- 
pressed in differential form (see also Binncy & Trcmainc 
1987). It follows that complete information about the ef- 
fect of isotropic softening on the equilibrium of 2-D discs is 
already contained in <S(|fc|s), which in this context has the 
meaning of a reduction factor for the transformed surface 
density. 

The relative quadratic deviations. The formulae for the 
angular speed fl(R) and k(R) are of particular interest 
since these quantities also determine the stability and re- 
laxation properties: 



n 2 (R) 



2irG 
~R~ 



H x 



k 2 (R) = 2ttG H 
4itG 



E(k)S(\k\s) (R) 



kt(k)S(\k\s) (R) 



(18) 



R 



Wi 



t(k)S(\k\s) (R), 



(19) 



S(fc) being an abbreviation for 7io[£(i?)](fc). The relative 
quadratic deviations of fl(R) and k(R) from their Newto- 
nian behaviours 



,(R) = 



\gHR)-gUR)\ 



(g = Q, k) 



(20) 



depend on the model of mass distribution and on s/Rd, 
i?d being the disc scale length. In order to estimate the 
magnitude of these deviations, we can set «S ~ 1 — \k\h 
and find that e g = 0(h/Rd) <C 1, which is strictly valid 
for R^$> h [the approximate formula for the circular speed 
is: v 2 (R) w v 2 N (R) + 2irGRT,'(R)h}. More important are 



the deviations that result near the centre. They imply a 
change in the number and/or location of the inner Lind- 
blad resonances, for a given pattern speed, together with 
larger natural scales for A, c, s, h and £rw- The contribu- 
tion of a massive bulge to the rotation curve makes the 
system more robust against such modifications, but in cer- 
tain respects this modelling is unclear and we are not yet 
in a position to draw quantitative conclusions (cf. Bertin 
1996; Bertin & Lin 1996; Junqueira & Combes 1996, Sect. 
2.2). 

2.2. 3-D models with isotropic/ anisotropic softening 

2.2.1. 3D vs. 2D 

The dynamics of 3-D discs with isotropic softened gravity 
is difficult to investigate because the effects of softening 
combine with those of vertical random motion in a compli- 
cated form. Nevertheless, we do in part understand how to 
benefit from the additional degree of freedom introduced 
into such models. A simple scale argument suggests that 
we should choose values of s sufficiently smaller than the 
Newtonian characteristic scale height, otherwise softening 
would significantly affect the vertical structure at equi- 
librium, i.e. both the mass distribution and the thickness 
scale. Likewise, the 2-D stability analysis suggests that in 
3D we should choose s <C s sa fe, otherwise softening would 
significantly interfere with the effect of thickness, as far as 
density waves are concerned. Note that choices of s > s cr ; t 
have the same consequences as in 2D, since the critical 
stability characteristics do not depend on the tempera- 
ture anisotropy. We mention that the vertical structure at 
equilibrium and the stability of 3-D discs with Newtonian 
gravity have been investigated in previous papers (Romeo 
1990, 1992; see also Paper I). 

The 3-D Jeans problem with isotropic softened gravity 
can be solved by using the technique of Fourier transforms 
(see, e.g., Sneddon 1972; Bracewell 1986). In Cartesian 
coordinates, this is the natural approach for factorizing 
the convolution in the Poisson equation and investigating 
the dynamical effects of softening, as Pfenniger & Friedli 
(1993) have previously emphasized. In the following dis- 
cussion we concentrate on the reduction factor and com- 
pare the Jeans problem with the dynamics of discs. The 
reduction factor for the volume density can be calculated 
analytically from: 



5j(|fc|a) = |fc|^ B [r^(r)](|fc|) ; 



(21) 



where k is the wavevector of the perturbation and J- s de- 
notes the Fourier sine transform: 



F s [g(r)}(\k\) = / g(r)sm(\k\r)dr. 
Jo 



(22) 



Equation (21) follows from the definition of Sj and the 
spherical symmetry of tp s , which allows expressing the 3-D 
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Fig. 2. The point- mass potential —Gmtp(r) (top) and the magnitude of the force Gm 2 \f\(r) (bottom) in isotropic (left) and 
anisotropic (right) types of softened gravity. In the isotropic types (discussed in Sect. 3.1) the gravitational interaction depends 
on the softening length s. In the anisotropic types (discussed in Sect. 3.2) the gravitational interaction depends on the planar and 
vertical softening semi-axes s\\ and s±, respectively. The anisotropic behaviour is only shown in the extreme cases (R, \z\ = 0) 
and (R — 0, \z\), and is labelled accordingly 



Fourier exponential transform in terms of JT S . Two alter- 
native equations that are more convenient for numerical 
computation can be obtained by first singling out the New- 
tonian behaviour of <p s (r) at long distances, and then inte- 
grating by parts (the integrand can be split in two ways) . 
As regards the magnitude of <Sj(|fc|s), Sj = 1 — o(|fe|s) at 
small \k\s and 5j <C 1 at large \k\s. The stability prop- 
erties are weakly affected by softening. In particular, sig- 
nificant deviations of the marginally stable wavelength \q 



from the Jeans wavelength Aj would occur for values of s 
comparable to Aj. For example, in the case of Plummcr 
softening this can be shown by using the simple numerical 
approximation Aq ~ Aj -I- |s for ^Aj ^ s < |Aj [the an- 
alytical formula for <Sj(|fe|s) is reported in Appendix A]. 
Weinberg (1993) has shown that in analogous models with 
Newtonian gravity the dominant source of noise is repre- 
sented by fluctuations on scales of stability interest. This 
fact and our stability analysis suggest that the collective 
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Fig. 3. The reduction factor S(k) of 2-D discs with isotropic (left) and anisotropic (right) types of softened gravity, k being 
the radial wavenumber of the perturbation. In addition, s is the softening length, s\\ is the planar softening semi-axis, and the 
types of softened gravity are abbreviated as in Fig. 2. Also shown are the ranges in which the contribution of self-gravity to the 
dispersion relation is stabilizing and destabilizing. (In the case with Newtonian gravity S — 1) 



relaxation properties are also weakly affected by soften- 
ing. Thus, there is a sharp contrast between the Jeans 
problem and the dynamics of discs: a flattened geometry 
strengthens the effects of softening, and rotation makes 
them critical for values of s that are an order of magni- 
tude below At- 



which allow expressing the 3-D Fourier exponential trans- 
form in terms of 7io and T c . Supposing that sii > sj_, then 
we expect that 5j(fc||,0) < <Sj(0, k±) for fcii = k±, which 
means that the stability and collective relaxation proper- 
ties are more affected parallel than perpendicular to the 
plane of reference. 



2.2.2. Anisotropy vs. isotropy 

The dynamical properties of 3-D discs are roughly decou- 
pled parallel and perpendicular to the plane, and a satis- 
factory modelling with isotropic softening may impose sig- 
nificantly different requirements on s. Pfenniger & Fricdli 
(1993) have discussed an ingenious way of coping with 
this difficulty, viz. to introduce a further degree of free- 
dom into such models: softening anisotropy. We assume 
that the point-mass potential —Grrvp(r) is of the form 
(f(r) — tp(R, \z\), and depends on the planar and vertical 
softening lengths sm and s±, respectively. The particular 
form in which softening anisotropy is implemented defines 
sii and sj_, and specifics their meaning in the context of 
the dynamical requirements. 

The 3-D Jeans problem shows the basic dynamical ef- 
fects of softening anisotropy. The reduction factor may 
be evaluated from (cf. 2-D discs with isotropic softened 
gravity): 



Sj(fc||,fej 



(*B 



fci)^c[£(fcn;|*|)](fc±), 



(23) 



where ip(ku;\z\) is an abbreviation for H.o[ip(R,\z\)](ku) 
and T c denotes the Fourier cosine transform. This equa- 
tion follows from the axial and planar symmetries of ip, 



3. Applications 

3.1. Isotropic types of softening 

As a first application of our method, we compare the dy- 
namical effects of three isotropic types of softened gravity: 
Plummer softening (Aarseth 1963; Miller 1970), cubic- 
spline softening (Hernquist & Katz 1989) and homo- 
geneous-sphere softening (Pfenniger & Friedli 1993); here- 
after abbreviated to P, CS and HS, respectively. Each 
name means that the softened point-mass potential can be 
viewed as the Newtonian potential generated by a spher- 
ical mass distribution of that type (but the force is eval- 
uated regarding the test particle as a point mass). In the 
following discussion we do not refer to this finite-sized 
particle interpretation of softening, unless otherwise spec- 
ified. The behaviours of these types of softened gravity are 
shown in Fig. 2 (left). In P the gravitational interaction 
is softened mainly at short distances and the Newtonian 
behaviour is recovered asymptotically, whereas in CS and 
HS softening is perfectly localized, i.e. beyond a certain 
cut-off distance the gravitational interaction is identically 
Newtonian. The form in which localization is implemented 
differs in the two types. In particular, CS has a rather soft 
cut-off at r = 2s, s being the nominal softening length, 
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Table 1. Stability and relaxation characteristics of 2-D discs with isotropic (top) and anisotropic (bottom) types of softened 
gravity 





Abbr. 


Softening Parameter 




Related Quantities 




Type of Softened Gravity 


a b 
*/opl ^?op 


'7op2 a 7/safc C 


»?crit 


C c 


op 


r) g 

Uop 


C/v h 


Acrit 


Plummer 
Cubic Spline 
Homogeneous Sphere 


P 
CS 

HS 


0.03 0.12 
0.06 0.23 
0.07 0.27 


0.25 0.20 
0.44 0.45 
0.52 0.53 


0.37 
0.61 

0.71 


1.70 
0.99 
0.83 


0.05 
0.11 
0.12 


0.39 
0.30 
0.30 


1.00 
0.44 
0.38 


0.37 
0.46 
0.47 


Homogeneous Oblate Spheroid J 
Homogeneous Prolate Spheroid J 


HOS 
HPS 


0.18 0.55 
0.03 0.11 


0.94 

0.24 0.18 


1.21 

0.37 


0.44 
1.73 


0.32 
0.05 


0.39 


0.04 
1.13 


0.61 
0.33 



a Values corresponding to one half the optimal relaxation level. 

b Optimal choice. 

c Safety threshold for approximate physical consistency. 

d Critical value. 

c Correction factor for the Rybicki- White relaxation time. 

f Optimal relaxation level. 

g Temperature anisotropy corresponding to the stability threshold for the optimal value of the effective thickness parameter. 

h Conversion factor, £ being the effective thickness parameter (below the safety threshold for approximate physical consis- 
tency). 

1 Critical radial wavelength. 

j HOS has s±/s\\ — i and HPS has s±/s\\ — 3, where s\\ and s± are the planar and vertical softening semi-axes, respectively, 
and sii is regarded as the softening length of reference. 



whereas HS has a sharper cut-off at r = s. The results 
of the dynamical comparison arc presented in Fig. 3 (left) 
and Table 1 (top), which are discussed below, and in Figs. 
4 (left) and 5. Useful analytical formulae are reported in 
Appendix A. 

Figure 3 (left) shows that the effective scale height h 
varies significantly with softening type, but apart from 
that the stability properties are analogous on scales larger 
than 2ft. This can be shown by rescaling k in terms of ft -1 
and noting that S is similar in the three types for \k\h < |. 
For \k\s <; 4, Shs changes sign and starts to oscillate, in 
contrast to Sp. An analogous behaviour occurs in SqSi 
but is less noticeable. Negative values of S mean that a 
given surface-density perturbation is in phase with the 
induced potential perturbation and, correspondingly, that 
the contribution of self-gravity to the dispersion relation 
is stabilizing. However, this feature does not affect the 
stability properties on scales comparable to the inverse of 
the typical radial wavenumber, since S changes sign well 
beyond the critical point /c cr its c rit [cf. Table 1 (top)]. The 
related oscillations of S mean that in certain ranges of \k\s 
noise is suppressed more effectively on larger than smaller 
scales, and correspond to oscillations of the short-wave 
branch of the dispersion relation. On the other hand, this 
feature has not significant consequences for the relaxation 
properties, because the corrected optimal relaxation level 
Cr op is approximately the same as in P [cf. Table 1 (top)] . 

Table 1 (top) shows that the stability and relaxation 
characteristics vary significantly with softening type. The 
largest variations occur between HS and P, concern Q/rj 



and 7/ sa fe, and are up to a factor of 3. In contrast, the cor- 
rected optimal relaxation level Cr op , S op and A cr it remain 
approximately constant. The variations between HS and 
CS are within 20%. 



3.2. Anisotropic types of softening 

As a second application of our method, we consider an 
anisotropic generalization of HS, the family of homo- 
geneous-ellipsoid softening (Pfenniger & Friedli 1993), 
and compare the dynamical effects of two representative 
spheroidal members: one oblate with sm : s± = 10 : 1, su 
and s±_ being the planar and vertical softening semi-axes, 
and the other prolate with s|| : sj_ — 1 : 3; hereafter ab- 
breviated to HOS and HPS, respectively. Regarding S|| as 
the softening length of reference, s_i_/ s || gives a measure 
of softening anisotropy. The behaviours of these types of 
softened gravity are shown in Fig. 2 (right). Anisotropy is 
implemented in a form consistent with the finite-sized par- 
ticle interpretation of softening, i.e. through a spheroidal 
deformation of the field particle. This mainly corresponds 
to a spheroidal transformation of the surface on which the 
force peaks, but it also strongly influences the behaviour 
in the plane. In particular, the degree of softening local- 
ization, the sharpness and magnitude of the force peak 
differ from those in HS [cf. Fig. 2 (left)], and are higher in 
HOS than HPS. So a 2-D analysis is required even though 
the isotropic case has already been investigated. The re- 
sults of the dynamical comparison are presented in Fig. 3 
(right) and Table 1 (bottom), which are discussed below, 
and in Fig. 4 (right). The major points are then general- 
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Fig. 4. The relative quadratic deviations en(R) (top) and e K (R) (bottom) of 2-D exponential discs with isotropic (left) and 
anisotropic (right) types of softened gravity, where £7 is the angular speed, n is the epicyclic frequency and the deviations 
are measured from their Newtonian behaviours. In addition, Rd is the disc scale length, s is the softening length, s« is the 
planar softening semi-axis, and the types of softened gravity are abbreviated as in Fig. 2. (In the case with Newtonian gravity 
en = e K = 0) 



ized in the final discussion. Useful analytical formulae are 
reported in Appendix A. 

Figure 3 (right) shows that, for fc|s|| Sa 2, <Shos is 
concave, in contrast to <Shps- The transition behaviour 
occurs in 5hs [cf. Fig. 3 (left)]. Concavity of S mainly 
means that softening has a stronger tendency to cause 
artificial stabilization for a given effective scale height h, 
for we know that the reduction factor of 3-D discs with 
Newtonian gravity is convex (cf. Paper I, Fig. 1). In fact, 
HOS does not mimic the effect of thickness for realistically 
large values of h, since the safety threshold for approxi- 
mate physical consistency s sa fc and the optimal value h op 
are ill-defined [cf . Table 1 (bottom)] . Concavity of S also 
means that softening causes a moderate degree of 'red- 
shift'. For \k\su ^ 4, 5hos an d ^hps become analogous to 



iShsj and analogous conclusions can be drawn concerning 
the stability and relaxation properties. 

Table 1 (bottom) shows that the stability and relax- 
ation characteristics differ considerably in the two types. 
The major qualitative differences concern ?7 sa f c and S op , 
and have been pointed out above [cf. discussion of Fig. 3 
(right)]. Excluding C/t], which has a restricted range of ap- 
plicability in HOS, the largest variation concerns r op and is 
beyond half an order of magnitude. In contrast, note the 
low sensitivity of the corrected optimal relaxation level 
Ct op , which remains roughly constant. 

Finally, instead of comparing those results with the 
isotropic case, let us clarify the dynamical relations be- 
tween all the spheroidal members of the family of homo- 
geneous-ellipsoid softening and the standard type of soft- 
ened gravity [cf. Table 1 (top)]. 
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Fig. 5. The reduction factor <Sj(|fc|s) for the 3-D Jeans prob- 
lem with isotropic types of softened gravity, where k is the 
wavevector of the perturbation and s is the softening length. 
In addition, the types of softened gravity are abbreviated as 
in Fig. 2. Also shown are the ranges in which the contribu- 
tion of self-gravity to the dispersion relation is stabilizing and 
destabilizing. (In the case with Newtonian gravity <S.j = 1) 



1. A general feature of the spheroidal members is that 
h is determined by s± alone: h = §s_i_. This is not 
intuitive since a more natural softening length in the 
plane is expected to be s\\. 

2. All the oblate members with s±/s« Ss I have an ill- 



defined s sa f c , i.e. s sa f c <; s cr ;t, and thus differ funda- 
mentally from P (cf. HOS). 
3. As s±/s\\ increases, such a difference becomes progres- 
sively less important (cf. HS and HPS). Indeed, the 
prolate member with s±/su = |, i.e. h = sn, turns 
out to be hardly distinguishable from P. (Values of 
s±/su *J 2-3 are not employed in simulations.) 

4. Conclusions 

Modelling gravity is a fundamental problem that must 
be tackled in A-body simulations of stellar systems, and 
satisfactory solutions require a deep understanding of the 
dynamical effects of softening. This problem has deserved 
special attention both in the past (e.g., Miller 1970, 1976) 
and in more recent times (e.g., Efstathiou et al. 1985; 
Hernquist & Katz 1989; Pfenniger & Fricdli 1993). Viewed 
in this general perspective, our contribution has a three- 
fold practical importance in addition to the points empha- 
sized in Sect. 1. 

1. The two present applications of our method reveal the 
dynamical differences between the most representative 
types of softened gravity: the isotropic P, CS and HS 



(abbreviated as in Sect. 3.1), and the anisotropic HOS 
and HPS (abbreviated as in Sect. 3.2). The major 
conclusions concerning their dynamical resolution, i.e. 
their faithfulness in simulating the Newtonian dynam- 
ics, are summarized below. 

(a) As regards the isotropic types, the dynamical res- 
olution is comparable. This results from the fact 
that, even though the spatial resolution and the 
effectiveness in reducing noise differ significantly 
in P, CS and HS for the same nominal softening 
length s, those differences can largely be removed 
by considering a more appropriate softening length 
of reference. 

(b) As regards the anisotropic types, the dynamical 
resolution is significantly coupled parallel and per- 
pendicular to the plane. In the plane, it decreases 
in quality from HPS to HOS, and the transition 
occurs in the oblate members for a softening ax- 
ial ratio s_l/s|| ~ § (HPS is dynamically similar to 
P). These disadvantages result from the finite-sized 
particle implementation of softening anisotropy. On 
the other hand, they have less importance than 
the advantage of introducing such a degree of free- 
dom into 3-D simulations of disc galaxies, which 
has been emphasized by Pfenniger & Friedli (1993) 
and in our method. 

(c) Last but not least, when employing these types of 
softened gravity in simulations of disc galaxies, we 
should recall that the dynamical resolution depends 
critically on two quantities: s, or s\\ for a given 
Sj_/sii, and the Safronov-Toomre parameter Q (cf. 
Paper I). The choice of s or s« should be checked 
vs. the profiles of the characteristic values s op , s sa f c 
and s cr it, which are tabulated in the applications. 
The choice of Q should be checked vs. the profiles 
of the stability threshold Q and level Q c j$, which 
can be evaluated as is explained in the method. 

Our method can be applied for testing new ideas about 
softening. There are two features that encourage such 
future applications. 

(a) One is the unified approach adopted for investi- 
gating stability, relaxation and equilibrium. As a 
result, full information about the dynamical effects 
of softening is contained in a single quantity: the 
reduction factor S. 

(b) The other is the modular structure of the method. 
We describe step by step how to extract detailed 
information concerning the dynamical properties, 
starting from S and pointing out the quantities of 
major interest. 

But our method can be applied in another, more fruit- 
ful, way: for developing new ideas about softening. In- 
deed, it opens a direct route to the discovery of optimal 
types of softened gravity for given dynamical require- 
ments, and thus to the accomplishment of a physically 
consistent modelling even in the presence of a cold in- 



Alessandro B. Romeo: Dynamical effects of softening in JV-body simulations of disc galaxies 



ff 



terstellar gaseous component. Such a future applica- 
tion will be the objective of a 'twin' paper. 
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A. Useful analytical formulae for Sect. 3 

S P (|fc| S )=c-l fe l s ; (Al) 

(CA?)p,CS,HS,HOS,HPS = lj to; 8> 80' 8 ' (^ 2 ) 

(%afe)p,CS,HS,HOS,HPS = 51 3T > 151 — 7 45 I (A3) 

(A4) 
(A5) 
(A6) 
(A7) 



> 70 ' 8 ' 80 ' 8 ' 

_ 1 14 8 8 . 

safeJP,CS,HS,HOS,HPS — 51 31 > 151 — , 45 , 

1 . 



(»7crit)p — e . 

(Acrit)p = - 5 

(n\ — 16 17 325 5 

(C)P,CS,HS ~ 3¥' T7504' 6 

5 J p(|fe|s) = |fe|sK 1 (|fe|s), 

sin 3 (|fc|s/2) 



5jcs(|fc|s) 



«5jHs(|fe|s) 



5 JHS (|fc|s/2). 



(|fc|*/2)3 

3 [sin(|fc|s) — \k\s cos(|fc|s)] 
(\k\sj* 



(A8) 



(A9) 



where K„ denotes the modified Bessel function of the sec- 
ond kind and order v (see, e.g., Abramowitz & Stegun 
1972). 

References 

Aarseth S.J., 1963, MNRAS 126, 223 

Abramowitz M., Stegun LA., 1972, Handbook of Mathemati- 
cal Functions - with Formulas, Graphs, and Mathematical 
Tables. Dover, New York 

Bender CM., Orszag S.A., 1978, Advanced Mathematical 
Methods for Scientists and Engineers. McGraw-Hill, Auck- 
land 

Bertin G., 1996, Spiral Structure in Galaxies: Competition 
and Cooperation of Gas and Stars. In: Block D.L., Green- 
berg J.M. (eds.) New Extragalactic Perspectives in the New 
South Africa. Kluwer, Dordrecht, p. 227 

Bertin G., Lin C.C., 1996, Spiral Structure in Galaxies - A 
Density Wave Theory. MIT Press, Cambridge 

Binney J., Tremaine S., 1987, Galactic Dynamics. Princeton 
University Press, Princeton 

Bracewell R.N., 1986, The Fourier Transform and Its Applica- 
tions. McGraw-Hill, New York 



Byrd G., 1995, Tidal Perturbations, Gravitational Amplifica- 
tion, and Galaxy Spiral Arms. In: Hunter J.H. Jr., Wil- 
son R.E. (eds.) Ann. N.Y. Acad. Sci. 773, Waves in Astro- 
physics. N.Y. Acad. Sci., N.Y., p. 302 

Combes F., Debbasch F., Friedli D., Pfenniger D., 1990, A&A 
233, 82 

Efstathiou G., Davis M., Frcnk C.S., White S.D.M., 1985, 
ApJS 57, 241 

Farouki R.T., Salpeter E.E., 1994, ApJ 427, 676 

Friedli D., 1994, private communication 

Fulbright M.S., Benz W., Davies M.B., 1995, ApJ 440, 254 

Gerber R.A., 1996, ApJ 466, 724 

Goodman J., Heggie D.C., Hut P., 1993, ApJ 415, 715 

Gradshteyn I.S., Ryzhik I.M., 1994, Table of Integrals, Series, 
and Products. Academic Press, Boston 

Gueron S., 1994, J. Comput. Phys. 110, 164 

Gurzadyan V.G., Pfenniger D. (Eds.), 1994, Ergodic Concepts 
in Stellar Dynamics. Springer- Verlag, Berlin 

Hernquist L., Barnes J.E., 1990, ApJ 349, 562 

Hernquist L., Katz N., 1989, ApJS 70, 419 

Hernquist L., Ostriker J. P., 1992, ApJ 386, 375 

Hockney R.W., Eastwood J.W., 1988, Computer Simulation 
Using Particles. Adam Hilger, Bristol 

Junqueira S., Combes F., 1996, A&A 312, 703 

Kandrup H.E., Smith H. Jr., Willmes D.E., 1992, ApJ 399, 627 

Masset F., Tagger M., 1996, A&A 307, 21 

Merritt D., 1996, AJ 111, 2462 

Miller R.H., 1970, J. Comput. Phys. 6, 449 

Miller R.H., 1976, J. Comput. Phys. 21, 400 

Palous J., Jungwiert B., Kopecky J., 1993, A&A 274, 189 

Pfenniger D., 1993, Order and Chaos in iV-Body Systems. In: 
Combes F., Athanassoula E. (eds.) iV-Body Problems and 
Gravitational Dynamics. Observatoire de Paris, Paris, p. 1 

Pfenniger D., Friedli D., 1993, A&A 270, 561 

Romeo A.B., 1990, Stability and Secular Heating of Galactic 
Discs. PhD thesis, SISSA, Trieste, Italy 

Romeo A.B., 1992, MNRAS 256, 307 

Romeo A.B., 1994, A&A 286, 799 (Paper I) 

Rybicki G.B., 1972, Relaxation Times in Strictly Disk Systems. 
In: Lecar M. (ed.) Proc. IAU Colloq. 10, Gravitational N- 
Body Problem. Reidel, Dordrecht, p. 22 

Schroeder M.C., Comins N.F., 1989, ApJ 346, 108 

Shapiro P.R., Martel H., Villumsen J.V., 1994, Adaptive 
Smoothed Particle Hydrodynamics and Galaxy Formation. 
In: Franco J., Lizano S., Aguilar L., Daltabuit E. (eds.) Nu- 
merical Simulations in Astrophysics. Cambridge University 
Press, Cambridge, p. 45 

Shapiro P.R., Martel H., Villumsen J.V., Owen J.M., 1996, 
ApJS 103, 269 

Shlosman I., Noguchi M., 1993, ApJ 414, 474 

Shu F.H., 1968, The Dynamics and Large Scale Structure 
of Spiral Galaxies. PhD thesis, Harvard University, Cam- 
bridge, Massachusetts, USA 

Sneddon I.N., 1972, The Use of Integral Transforms. McGraw- 
Hill, New York 

Sommer-Larsen J., Vedel H., Hellsten U., 1997, ApJ (submit- 
ted) 

Theis Ch., 1997, in preparation 

Toomre A., 1963, ApJ 138, 385 

Vandervoort P.O., 1970, ApJ 161, 87 



12 Alessandro B. Romeo: Dynamical effects of softening in iV-body simulations of disc galaxies 

van Veldhuizen M., Nieuwenhuizen R., Zijl W., 1994, J. Corn- 
put. Phys. 110, 196 
Weinberg M.D., 1993, ApJ 410, 543 
Weinberg M.D., 1996, ApJ 470, 715 
White R.L., 1988, ApJ 330, 26 
Zhang X., 1996, ApJ 457, 125 



This article was processed by the author using Springer- Verlag 
DTeX A&A style file L-AA version 3. 



0) 
00 





- 


o 


1 1 1 1 1 

O 


\ 1 1 1 . 

\ .-■' 
\ J 
\ / 






\ 


\ 


\ /- 

\ / 






ST 


\Qf 


\ .■ / 

A / 

\ / 


Ml 


- 






\ / 
\ / 


fl 








/ 












n 


" 






/ \ 


<u 








/ \ 


^> 










*fH 










o 








. 


CO 








\ 


K 


- 






\ 


H 








\ 


s 


~ 






\ 


S 








\ 


£> 








\ 


j 








* 


cu 


- 




, , , 1 , 


\ - 

\ 

\ 



LO 

o 
d 



d 



d 



0) 
00 











' 




1 




; 








^^B9 




M 






Kti 


Eek3 








- 


c 






V*J 






























c 




1^1 














II) 
















- 


4J 


















e^ 

























H<H 


kriH4 










Ul 




Brl 




Es^Efl 


El 








«fH 






Rifl 


^H?il 


K*S 








o 




■SI 


^El 


Hr^^lsl 


eeI 




F ' 


- 


CO 

0) 




■a 


HE] 


H*^H£fl 


Ejfi 




/ ■'' 


_ 


ft 














/' ■' 




>, 














/// ■' 


- 










H 






GO 


Ul 






/// ■ 




LJ 






o 


Oh 






/// ■' 




Dh 




^ 


K 


K 






:'// ■■' 


- 


< > 




i 


i 








// / 




K 






i 






£ //; 


■ 


H 
O 




1 


i 






^^/ 




- 


Ul 










,' 


^ — / 


■ o 




< 




i 


/,-- 

u 


^ 


''' 


o" 

1 i 


i i i 


- 



_N_ 

a; 











| , , , 


| i i i i ; 








Kl 
O 


CO 


f- 






z 


ffi 


K 


I: _ 


Q0 







































a 










N ■ 




<u 












I. - 


+j 










O 




4^ 














o 














ui 














^ 












' «: 


o 












'■/' : - 


Ul 










'/' ' 


<u 










ft : — 


Ph 










//:■' 


X 










//''•* " 


E-h 












U 










'X ' ■' 














(X 

o 








o 


/'/ ' ■'' 


K 








S 


---'"' '^/r / : ~ 


E^ 










^^^^ ''' - — 


O 


*." 


OT 








^^^.^^ 


"^^ ,''' o" \ U 


<; 








— """"--.^ 


,'"' S ''■■ 


















~n~ 




__ — •"* 


""* -* -. 


















o 




„ — — "" 


"- -^ 

































"-, 








0) 
00 





1 1 1 1 1 1 1 1 1 


1 
\ 
\ 






o 




. 




m go 


g 




\ 




- O P-, 




\ 


Ofi 


ffi K 


E3 




I 


ti 


— 1 : 


3 




\ _ 




: 


S 






CI 


: 


< 






0) 


1 : 


E-h 






^ 




m 






<+H 




K 




I 


o 




Q 




/ 


in 








/ 


o 




H-> 

> 




/ 
/ 


m 




cd 


/ 








dJ 






' 




O 


— 


a 






^M 


g 








'qj 


: 1 


E3 

3 


- 






go 








U 








S 




E 








O 




en 


_ 


K 


,'' 






E- 


,' 


>> 




O 
W 


s' 


> 

CO 




s 


,''' 


o 


- 


< 


,' 








— / 




— 






'a; 






/ ..-■'■ 


CO 






' 1 
...■■1 1 1 1 1 1 1 1 1 1 1 



o 



m 



^ 



•a. 



cv 



^ 





l 1 l l | 1 l 1 l 








a 










g 








go go 










Ph u K 


R 






&0 


— 1 1 


3 






a 


1 


5 










<* 






'S 


1 


h 






0) 




CO 






H-> 




H 






H-H 




Q 






o 










w 




>> 




: 


o 




> 

ca 








~ 


Sh 

o 




g 




ft 




»^ 


■i 


M 








1j 
GO 


■•''/ 


3 
5 
<; 




U 


' ' ■' ' 1 








h 




E 


.■■' / / 


CO 


_ 


o 


■' '' / 






K 


■ ' / / 


H-3 




E- 


■■' '' / 


> 
CO 




O 


..-■' ^' / 




CO 


.-'' s' / 


Eh 






■'-''' // 


C3 
























CD 








GO 













































d 



"* 



•a. 



^ 



0) 
00 







1 ' 1 


1 ' 


j 




m 


m 




i 




- o 


P-, 




j 


on 


X 


K 




: 1 


C 


- 






■ 1 
: 1 


C 












_ 






i 


«fH 










<) 










m 








41 












o 










0) 








1 


ft 










>-, 










H 










U 








T 


Oh 










C ) 








— [■ 


K 










H 










O 










in 


*"H 






J 


s 


o 






1_ 


< 


a; 
\ 


1 ..-,--1 


y .. 


1 



a; 



LO 

d 




Q. 
Q, 







' 


1 






m 


GO 








- o 


Oh 




_ 


ao 


ffi 


X 






a 


- 






; "I 


"a 










CD 








-1 


^-> 










<+H 










o 








1 


w 










<fH 










o 


_ 






" 


w 










CD 


- 






-| 


a 








1 


>^ 








1 


E-i 






1 
1 


U 






1 

-f 










PL, 








o 


. 




_4. 


K 








O 


- 




1 


W 


*"H 






s 


_ d 

1 




f 


< 


Tj 








a; 




r 




vr 


1 I.------I' 


-lu'— J- — T — T"'i 



a; 



in 
d 



CO 

d 



CD 

d 



d 



d 



m 





- 








t 

i: 
f 






CO 


w 




I: 




Ph 


o 


ffi 




I 1 


ao 


_ 


1 






1; 


Pi 




1 






i; 














Pi 


_ 


1 






i- 


0! 










i: 


+J 












«fH 












o 










I: 


w 










i; 


M-) 


- 








1: ~ 


o 










i; 


w 


- 








I: ~~ 


<u 










i: 


ft 










';■ 


>» 












E-i 










: 


U 








/ 
/ 












/ i 




O, 









/ i 





O 








/ i 




K 








/ '■ 




E- 


" 






/ '■' 


~ 


O 


. 






/ '•' 




W 














d 

II 






/ '•'' 


- 




d 




-j=-^r 


^^r«°f'\ 








d 



CO 

d 



CD 

d 



d 



d 



0) 
00 





i 




| i i i i 











o 












g 








CO 


Ul 










Ph u 


X 


R 




■ 


&0 


I ! 




3 




■. 


Pi 


| 




5 












<* 






fl 


_ | 




h 






0) 






CO 






H-> 






H 






^H 






Q 




o 








/' 


w 






>> 


/' ■■' 


o 






> 

ca 












.'// 


g 




ft 






1 


.-'■ / 


M 












■■■' / / 


m 

<! 




O 






.-•■"' / / 












h 




a, 






.,-■' s' / 


OT 


_ 


o 






s / 






K 






s' / 






o 




■"' /'' 




> 




CO 








!h 






.■' y' 






1 




































cu 






// ~^~ 






m 










1 I 



"* 



_w_ 

.42 



d 



